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Topological physics in photonic systems have attracted great attentions in recent years. In this work, we theoret¬ 
ically study the one and two dimensional photonic quasicrystal resonator lattices characterized by the first and 
second Chern number, which show exotic boundary states within the photonic energy band gap. In particular, 
the second Chern number protected edge states has opened up new possibilities for realizing topological physics 
of dimensions higher than three in photonic systems, which is highly sought for. Such photonic systems can be 
easily experimentally realized in regular photonic crystal with dielectric rods in air, by varying the radius of 
the rods, so we propose experiments realizing our predictions. 
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1 Introduction 

Topological phases in condensed matter including quantum Hall effect (QHE)^, topological insulators^”^, and quan¬ 
tum anomalous Hall effect (QAHE)^“^ have recently garnered tremendous interest in condensed-matter physics, ma¬ 
terial science and electrical engineering communities^^. In those systems, each energy gap is attributed an topological 
index. A nontrivial topological index is usually associated with interesting boundary states, quantized response, and 
exotic quasiparticle excitations. Soon after the discovery of these topological phases in condensed matter, similar topo¬ 
logical effects in photonic systems have been proposed and demonstrated, in analogy to the topological physics 
in condensed matter. Although topological phases exist in any dimension, those systems have dimensions below 3D, 
the maximum real space dimension in our universe. Thus for topological phases described by Chern number, only 
2D systems protected by the first Chern number including QHE, QAHE and its photonic analogy are realized. 
Topological phases protected by the second Chern number exist in 4D^^, and thus were thought to exist in theory only. 

Only recently, using the physical properties of quasicrystals (QCs)- nonperiodic structures with long-range 
order, 4D QHE protected by second Chern number have been proposed in 2D photonic systems^^. This is due to QCs 
can oftentimes be derived from periodic models of higher dimensions, so they contain extra information in higher 
dimensions^^”^^. The essential physics in their proposal is the on-site modulation of energy However, in this 

proposal, z dimension is treated as time and thus the wave vector in z, (3z represents energy. This is different than the 
conventional photonic crystal (PC) context, in which the frequency uj is the energy and its experimental setup has been 
very mature to realize exotic photonic physics. 
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It is known that the resonant frequency of resonators in photonic crystals can be modulated through their size^^. 
With this principle, in this work we propose a 2D photonic crystal resonator lattice that have second Chern number 
protected edge states within the band gap of frequency uj . The proposed structure can be easily manufactured with 
modern fabrication technology. 


2 Results 


We first introduce some basic theory on the tight binding formalism of resonators in PCs. With that, we first illustrate 
our idea with first Chern number protected one dimensional PC resonator boundary states, then present our main 
results on the second Chern number protected two dimensional PC resonator boundary states. 


Resonator Defect States in Photonic Crystal In this work, we consider the square lattice PC composed of dielectric 
rods with radius 0.20a and refractive index 3.4 in air (Fig. Its TM modes are well known to have a band gap 

from 0.29 — 0.42(27rc/a)^^. Like dopant in a semiconductor, when the radius of a rod varies from 0.20a, the defect 
states form localized states in the photonic band gap. Using a hamiltonian formalism^^, where the resonator mode 

is the eigenvector and the resonant frequency uj is the eigenvalue of the equation, we have 

Qo lip) = UJ lip) ( 1 ) 



with 


00 


0 isQ^Vx 

—i/ig^Vx 0 


( 2 ) 


Here uj is the resonant frequency lying within the bulk PC band gap (Fig. [^). So{x, y) is a function of the refractive 
index in real space. Inside the dielectric rods 5o = 3.4, outside = 1. When the radius of the defect r < 0.2a, the 
defect states are monopole states (Fig[^). 
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Figure 1: (Color online) (a) The one dimensional resonator array (small white rods) embedded in two dimensional 
photonic crystal of dielectric rods (big white rods) with radius 0.20a and refractive index 3.4 (AlGaAs). The resonator 
defect is created by changing the radius of a single rod. For example, the case where r = 0 corresponds to the removal 
of a rod. The green box indicates the unit cell we use to simulation the energy band diagram of the whole structure, 
including both the 2D photonic crystal and its embedded defect array. (b)The Brillouin zone of the unit cell in (a), (c) 
The energy band diagram of the whole unit cell in (a) in TM mode, the dielectric rods (big white rods) with radius 
0.20a form a band structure with band gap 0.29 — 0.42(27rc/a), inside the band gap are the tight binding states of the 
ID resonator defect array with radius 0.0965a.(d) The electric field distribution of the defect states in TM mode. The 
defect states are monopole states and form a tight binding lattice. 


One Dimensional Photonic Crystal Resonator Lattice When we have more than one resonators that are close 
enough in space, they are no longer isolated and the perturbed Hamiltonian is now 0 = ©o + 4^, where 


0 i{e ^ — Sq ^)Vx 
0 0 


The perturbation terms result in 




- / AeuJoE^Ei^dV 
2^ ^eo\E^\^dV^ ^eo\E^\^dV 


(3) 

(4) 
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where = e — £o is the difference of spatial dependent refractive index of the unperturbed system (single defect) eo 
and the perturbed system (coupled defects) Vaa is the self-energy correction for ccq of defect and Vai 3 is the 
hopping between defect a and Defining cCq, = ccq + Vaa, we can write the hamiltonian as 

H = E- a iV^a) I ^ ^ Va(3 iV^a) ("0/31 (5) 

ot a,13 

For an ID resonator lattice of lattice constant 3a with translational symmetry (Fig. considering nearest neighbor¬ 
hood hopping only, Va,a-\-i = Va-i,a = t,uja = ^c- Using Bloch theorem, we can solve the energy to be 

uj{k) = ccc + 2tcos{3ka) (6) 

as shown in Fig. [^. By fitting Eq. we can get for r = 0.0965a, ujc = 0.3312231(27rc/a), t = —0.00485(27rc/a). 
Fig. shows the Wc for lattice of resonators when we change the radius. 

One Dimensional Photonic Quasicrystal Resonator Lattice We have discussed ID photonic crystal resonator lat¬ 
tice composed of resonators with translational symmetry. Next we will elaborate how to construct a quasicrystal 
resonator lattice with bulk bands characterized by the first Chem number. Here the word “quasicrystal” does not mean 
a variation of lattice spacing, but a variation of each resonator’s radius. The variation of resonator’s size introduces the 
on-site modulation of energy needed for nontrivial Chern number^^“^^. Here are the construction steps: 

(1) Using ujco of r = 0.0965a as the zero energy reference, we compute the on-site energy of resonator site x: 

= ^co + cos(27r5x + 5 is an irrational number. In this paper, we use 6 = (1 -1- \/5)/2. x is an integer index 

for sites. 

(2) Using the a monotonic function of ujx of the resonator radius in ID resonator lattice with translational 
symmetry discussed in last section as shown in Fig. |^, we can determine the radius Vx for resonator x by inverting 
the function. 

Now we get a resonator QC lattice, each site with a different size. Notice that because ujx in step (2) was 
computed assuming ID translational symmetry, now for the quasicrystal lattice without translational symmetry, uJx is 
modulated by its neighborhood with different sizes, same for the hopping parameter t. Considering nearest neighbor¬ 
hood hopping only, now the hamiltonian is 

^ = E- X \'^x) {'^x \ ^ ^ Vx^x+l \'^x) (V^cc+ll “1“ H.c. (7) 

X x^x-\-3 

whereas we have defined functions uJx = uj{{27Tbx -h (l)x)raod27r), Vx^x+i = Vi{{27Tbx + (j)x)mod2TT), Vx^x-i = 
V2{{27Tbx + (t)x)mod27T) . These functions can be defined because according to our steps to construct the lattice, 
{27Tbx + (j)x)mod2T: determines the size of the resonator at x as well as its neighborhood. Next we shall proceed to 
show Eq. |7] describe a Chern insulator with first Chern number. 
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In order for the Chern number to be defined, we introduce a twisted boundary condition by parameter 

H = y] \'>Px) {'>Px+l\+ X] IV’a:) (V’a:-l| ( 8 ) 

X x^x-\-l x,x—l 

Next we set 5^ = [6 • LJ / L as a rational approximation to b so that the modulation is periodic. When L +oc 
the results will be the same as irrational b. When L is a prime number, for a shift (p ^ (p e where 5 = 21^1 jL 
and / = 0,1,L — 1, we can always find a translation of lattice site n ^ n ^ GO, 1,L — 1, so that any 

function V{{27rbLn + 0 + s)mod27T) = V{{27rbL{n + n^) + (j))mod2T:). We denote the translation operator by 
sites as T^, then + 5 ) = T^Hi){(j))T~^. The equivalence between phase shift and translation means that in the 
thermodynamical limit L + 00 , the shift has no physical consequence such as gap closing and the energy spectrum 
is independent of To evaluate the Chem number we define the projector operator for states below the band gap as^^ 




X] |n) (n| 


Er^<E 


gap 


(9) 


\n) is the state with energy and Egap is any energy within the band gap. Since the projection operator differs 
from the hamiltonian only by its eigenvalues, we also have P(0 + 5 , 6>) = TsP{(l)^0)T~^ and d(f)P{(l) e^O) = 
T^d(f)P{(j)^ 0)T~^. The Chern number can be calculated as 


Ci = — [ d<t)dec{e, (f>) 
27^^ Jo 

where 


C{cP,e)=Tr[P 


dP dP 
Ei>'~d9 


is the Chem number density^®. It can be easily proven that 


( 10 ) 


( 11 ) 


c{(j) + e,e) = c{<j>,e) 


( 12 ) 


So C(0, 0) has 1/L periodicity. Thus 

p 27 t p 27 tIL 

/ 0) = L d(t)C{(t), 0) = 27rC(0 = 0,6>) + 0{l/L) (13) 

Jo Jo 

Next we need to prove C(0, 0) is also independent of 0. From Eq. « 1 P iOjL. According to perturbation 

theory, H and P are suppressed by a factor of l/L, thus in the thermal dynamic limit, C(0, 0) is independent of 0 as 
well, and 

Cl = -{27ri)C{(t) = 0,6> = 0) + 0(1/L) (14) 

Eq . |7]can be well approximated as (See Methods Section) 

H = {UcO + K COs{2TTbx + (j):c)) l-fpx) {i’xl + E * IV’®) (V’a;+ 1 1 + E (^5) 

x,x-\-l x,x—l 

The bands of this model with periodical and open boundary condition are plotted in Eig. ^ and c, with the number of 
sites chosen to be 68, t = —0.00485 and A = 0.010412. Here according to Eig. A is determined by the range of 
resonator radius in the QC lattice. Eor the value we use, Vx varies from 0.0875a to 0.1065a. By comparing those two 
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figures, we can clearly observe the boundary states appearing in the open boundary condition, whose number equals to 
the Chern number difference of the bulks bands across the gap. In Fig. we plot the wave function of the boundary 
state at = 1.174 and = 0.3361(27rc/a) using the tight binding basis at different sites x, which can be simply 
experimentally realized using the construction steps described in this Section. 





Figure 2: (Color online) (a)The energy at F point (Fig. [^) of the ID resonator defect lattice of various sizes embedded 
in an 2D PC of dielectric rods with radius 0.20a (Fig. [T^). (b) The bands of the ID defect QC (Fig. [T^)constructed 
in Sectionj^with periodical boundary condition according to Eq. the number of sites are chosen to be 68, and the 
parameters are t = —0.00485, A = 0.010412, 6 = (1 + \/5)/2.(c) The bands of the same ID resonator QC as in (b) but 
with open boundary condition instead, (d) From (c), we can observe that when (j)x = 1.174 and = 0.3361 (27rc/a), 
only boundary state exists. Here we plot the wave function of this state using the tight binding basis at different sites 

X. 


Two Dimensional Photonic Crystal Resonator Lattice The hamiltonian to describe two dimensional photonic crys¬ 
tal resonator lattice is still Eq. whereas a = (x, y). For a 2D resonator lattice of lattice constant 3a with transla- 
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tional symmetry, considering nearest neighborhood hopping only, V(^x,y),{x+i,y) = ^(cc,y),(cc-i,y) = y{x,y),{x,y-i) = 
y{x,y),{x,y+i) = ^{x,y) = ^c- Using Bloch theorem, we can solve the energy to be 

uj(k) = cjc + 2t{cos{3kxa) + cos{3kya)) (16) 


as shown in Fig. [^. By fitting Eq. 16 we can get for r = 0.095a, cjc = 0.3280711 (27rc/a), t 
Fig# shows the Wc for lattice of resonators with different radius. 


-0.0051005(27rc/a). 



Figure 3: (Color online) (a) The two dimensional resonator array (small white rods) embedded in two dimensional 
photonic crystal of dielectric rods (big white rods) with radius 0.20a and refractive index 3.4 (AlGaAs). The resonator 
defect is created by changing the radius of a single rod. For example, the case where r = 0 corresponds to the removal 
of a rod. The green box indicates the unit cell we use to simulation the energy band diagram of the whole structure, 
including both the 2D photonic crystal and its embedded defect array. (b)The Brillouin zone of the unit cell in (a), (c) 
The energy band diagram of the whole unit cell in (a) in TM mode, the dielectric rods (big white rods) with radius 
0.20a form a band structure with band gap 0.28 — 0.42(27rc/a), inside the band gap are the tight binding states of 
the 2D resonator defect array with radius 0.095a. (d) The electric field distribution of defect states in TM mode. The 
defect states are monopole states and form a tight binding lattice. 
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Two Dimensional Photonic Quasicrystal Resonator Lattice We have discussed 2D photonic crystal resonator lat¬ 
tice composed of resonators with translational symmetry. Next we will present the main results of our paper on how to 
construct a QC resonator lattice with bulk bands characterized by the second Chern number in a way similar to the ID 
QC case. The variation of resonators size introduces the on-site modulation of energy needed for a nontrivial Chern 
number^^. Here are the construction steps: 

(1) Using ujco of r = 0.095a as the zero energy reference, we compute the on-site energy of resonator site 

(x, y): = ^co + cos(27r5x -h (p^) + \ cos(27r5^ -1- py). b is an irrational number. In this paper, we use 

6= (1 +V5)/2. X is an integer index for sites. 

(2) Using the a monotonic function of uj(x,y) of the resonator radius in the 2D resonator lattice with translational 
symmetry discussed in last section as shown in Fig. |^, we can determine the radius r(^x,y) for resonator (x, y) by 
inverting the function. 

Now we get a resonator QC lattice, each site with a different size. Notice that because uj(^x,y) in step (2) was 
computed assuming 2D translational symmetry, now for the QC lattice without translational symmetry, uJ(x,y) is mod¬ 
ulated by its neighborhood with different sizes, same for the hopping parameter t. Considering nearest neighborhood 
hopping only, now the hamiltonian is 

H ^{x,y) \p{x,y)) {'P{x,y) \ + ^{x,y),(x±l,y) \'P{x,y)) {'P{x±l,y) | + ^{x,y),{x,y±l) \'P{x,y)) {'P{x,y±l) \ 

x,y x,y x,y 

(17) 

whereas we have defined functions = oj{{2'Kbx+4>x)mod2T:, {2'!Tby+(py)mod2Tr), V(^x,y),(x +i,y) = Ui((27r5x+ 

px) mod27T, {27Tby + py)mod27r), = V2{{27rbx + px)mod27T, {27rby + py)mod27r), V(^,^),(a:,y+i) = 

Vs{{27rbx + 0a^)mo(i27r, {27rby + py)mod27T), V(^x,y),{x,y-i) = V4((27r6x + px)'^od27T^ {27Tby + py)mod27r). These 
functions can be defined because according to our steps to construct the lattice, {27Tbx + 0a^)mod27r and {27rby -1- 

py) mod2TT ubiquitously determine the size of the resonator at {x^y) as well as its neighborhood. Next we will 
elaborate how Eq. [^describes a Chern insulator with second Chern number following the ID QC case. 

In order for the Chern number to be defined, we introduce a twisted boundary condition by parameter 0^^ 

H = ^^^(x,y) \'P(x,y)} {p(x,y)\ ~^'^^'^{x,y),(x-\-l,y)^ ^ \'P{x,y)) {P{x+l,y) \ + ^{x,y),{x-l,y)^ \P{x,y)) {dP{x 

x,y x,y x,y 

+ '^^^(x,y),(x,y-\-l)^ \P{x,y)) {dP{x,y+l)\ ^'^^^^y{x,y),{x,y-l)^ \P{x,y)) {dP{x,y-l)\ 

x,y x,y 

In the same procedure as in the ID QC case, by setting = [5 • TJ/L as a rational approximation to b we can 
prove the equivalence between phase shift and translation in the thermodynamical limit L + 00 . So the shift has no 
physical consequence such as gap closing and the energy spectrum is independent of p. To evaluate the Chern number 
we define the projector operator for states below the band gap as 

P{(j)x,<py,0x,0y) = ^ \n) {n\ (19) 


I 

( 18 ) 
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\n) is the state with energy and Egap 
Chern number can be calculated as 


is any energy within the band gap. Given 


^2 = 


p27T 

Jo ' 




where 


al3^6 


dP d ^dP dP\ 

- ——i r I Jr' - P - 

\ d(t)a 94)13 94)^ 94)5) 




is the Chem number density. As in the ID quasicrystal case, 


{4>:c, 4>v,dx,0y), the second 

( 20 ) 

( 21 ) 


c{4 + s,0) = c{4),e) 


( 22 ) 


and C{4), 0) is also independent of 9 in the thermal dynamic limit. 

C 2 = (27r)4C(0^) ^ 0 


(23) 


Eq. 17 can be well approximated as (See Method Section) 


H = ^(wco + Aa;Cos(27r6a; + ?)^) + Ay cos(27r&y+ (/)j^)) \'4){x,y)) {'>P{x,y) \ l^(^-y)) {'^(x+i,y)\ 

x,y x,y 

“1“ ^ V ^ ) (,'^(x,y-\-i) I E H.c. (24) 

x,y 

The bands of this model with periodical and open boundary condition are plotted in Fig. and c, with the number of 
sites in both x and y directions chosen to be 34, t = —0.0051 and Xx = ^y = 0.00741925. For the value of A we use, 
r(^x,y) varies from 0.081a to 0.108a. By comparing those two figures, we can clearly observe the edge states appearing 
in the open boundary condition. In Fig. we plot the wave function of the boundary state at (jy^ = —2.443, = 0 

and ujx = 0.3427(27rc/a) using the tight binding basis at different sites (x, y), which can be simply experimentally 
realized using the construction steps described in this Section. 
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Figure 4: (Color online) (a)The energy at F point (Fig. [^) of the 2D resonator lattice of various sizes embeded in 
an 2D PC of dielectric rods with radius 0.20a (Fig. |^). (b) The bands of the 2D defect QC (Fig. [^)constructed in 
Section with periodical boundary condition according to Eq. the number of sites in both x and y directions are 
chosen to be 34, and the parameters are tx = ty = —0.0051, Xx = = 0.00741925, 5 = (1 + \/5)/2, = 0.(c) 

The bands of the 2D defect quasicrystal the same as in (b) but with open boundary condition, (d) From (c), we can 
observe that when = —2.443, = 0 and ujx = 0.3427(27rc/a), only boundary state exists. Here we plot the wave 

function of this state using the tight binding basis at different sites (x^y). The boundary state is localized on one edge. 


3 Methods 

Topological Equivalence between ID Models First we find that the hopping term in Eq. [^is a special case of Eq. 
between sites i, j, so we just need to consider the generic case, defining Afij = — ro, , where 

vq is the average of the resonator radius 0.0965a. Also we use first order approximation — vq = c{ujx — cCco) = 
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cXx cos{2TTbx + (f>x), whereas c is a constant 


Vi 


j = V(ri,rj) 

= V{Afij,Arij) 

dV d'^V o 

dV d'^V o 

dV cA 


d^V 


V^(0,0) + ^(cos(27r&i + (p^) + cos(27r&j + (p^)) + 


cA 

- ^ aAf- + cos(27r6j + 0^)) 

= t + Fi(cos(27r6i + + cos(27r6jf + 0cc)) 


(25) 


The first order term of vanishes because according to Eq. Vij is invariant under the exchange of i,j. In the 
second last step we ignored the second order term. Similarly 


= uj{rx, 0, 0) + -Ara,,a:-1 + 


duj 


duj 


-^^x,x-l + C- 


^Ar x^x+i 
duj 


Ar. 


'x,x-\-l 


-Acc. 


x,x-\-l 


dArx,x-i ’ dAvx^x+i 

= CCa: + F2{cOs{27rbx + — COs(27r6(x — 1) + 0a^)) 

+ F2(cos(27r6x + 0^^) — cos(27r6(x + 1) + 0cc)) 


(26) 


According to numerical simulation, Fi and F 2 are negative and positive and bounded by Fm = 0.0001492625 <C 
|t| , I A| in Eq. 15 Thus Eq. [^can be very well approximated by Eq. 15 We can also see their topological equivalence 
by plotting the energy states at = 0 as a function of F = Fi = — F 2 = —aFm, where a G [0,1] as a control 
parameter (Eig. Because the bulk band gaps did not close as F evolves from 0 to —F^, the two models are 
topologically equivalent. 
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Figure 5: (Color online) (a)The evolution of bulk bands of Eq with first order approximation as a parameter of F 
defined in Eq. |^and|^(b) The evolution of bulk bands of Eq with first order approximation as a parameter of F 
defined in Eq. |27|and|^ 


Topological Equivalence between 2D Models Same as the ID case, the hopping term in Eq. [T^is a special case of 
Eq. [^between sites i, j, so we just need to consider the generic case, defining — ro, Ar^^j = 

where ro is the average of the resonator radius 0.095a. Also we use first order approximation rx,y — ro = c{uj(^x,y) ~ 
<^co) = cos(27r6x + (px) + ^y cos{27rby + 0^)), where c is a constant 


= 


dV d‘^V 

dV d^V 

dV c, „ , 


(27) 
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Similarly 


^ix,y) — ^{'^{x,y)^ ^'^{x,y),{x — l,y)j ^'^{x,y),{x-\-l,y) j ^'^{x,y) ,{x,y—l) ^ ^'^{x,y),{x,y-\-l)) 

duj duj 

= ^{'^{x,y)^^^ 0: 0: 0) + ^ A ~ -^^{x,y),{x-l,y) + ^ a - ^^{x,y),{x-\-l,y) 

(at', .in A T, — 1 . 7/1 Lfi-Xi (■ 


dAr 


^(x,y),(x-l,y) 

duj duj 


{x,y),{x,y-l) 


dAr 


{x,y),{x,y+l) 


\x,y),{x+l,y) 
^'^{x,y),{x,y+l) 


duj 


^(x,y) ~l~ ^ 
+ C 


9^'^(x,y),(x-l,y) 
duo 


Auj 


duj 


{x,y),{x-l,y) 


SAr, 


-Auj. 


{x,y),{x+l,y) 


{x,y),{x+l,y) 




duj 


-Auj 


{x,y),{x,y+l) 


{x,y),{x,y+l) 


(28) 


With Aa; = \y, we can define Fi and F 2 the same as the ID case. According to numerical simulation, Fi and F 2 


are negative and positive and bounded by = 0.00019138 ^ |t|, |A| in Eq. 24 Thus Eq. 17 can be very well 
approximated by Eq. 24 We can also see their topological equivalence by plotting the energy states at = 0 

as a function of F = Fi = —F 2 = —aF^, where a G [0,1] as a control parameter (Eig. [^). Because the bulk band 
gaps did not close as F evolves from 0 to —F^, the two models are topologically equivalent. 
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